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Abstract 

In this paper, a multiscale modelling strategy is 
used to study the effect of grain-boundary sliding on 
stress localization in a polycrystalline microstructure 
with an uneven distribution of grain size. The 
development of the molecular dynamics (MD) analysis 
used to interrogate idealized grain microstructures with 
various types of grain boundaries and the multiscale 
modelling strategies for modelling large systems of 
grains is discussed. Both molecular-dynamics and finite- 
element (FE) simulations for idealized polycrystalline 
models of identical geometry are presented with the 
purpose of demonstrating the effectiveness of the adapted 
finite-element method using cohesive zone models to 
reproduce grain-boundary sliding and its effect on the 
stress distribution in a polycrystalline metal. The yield 
properties of the grain-boundary interface, used in the FE 
simulations, are extracted from a MD simulation on a 
bicrystal. The models allow for the study of the load 
transfer between adjacent grains of very different size 
through grain-boundary sliding during deformation. A 
large-scale FE simulation of 100 grains of a typical 
microstructure is then presented to reveal that the stress 
distribution due to grain-boundary sliding during 
uniform tensile strain can lead to stress localization of 
two to three times the background stress, thus 
suggesting a significant effect on the failure properties of 
the metal. 

Introduction 

Uneven stress distribution and stress localization 
during deformation are the key factors for fracture and 
failure in polycrystalline metals. The inhomogeneous 
polycrystalline microstructure that consists of grains of 
different size and shape joined together at different angles 
and forming various types of grain-boundaries (GBs) 
creates inhomogeneous deformation fields at, otherwise, 


homogeneous loading. There are a number of factors 
responsible for the appearance of inhomogeneous 
deformation inside the polycrystal. The co-existence of 
grains of different size and orientation with anisotropic 
elastic properties is one factor. The difference in the 
structure and properties of the GBs between grains of 
various missorientations is a second factor for uneven 
stress distribution. A third factor is grain-boundary 
sliding (GBS), i.e., the rigid translation of one grain 
relative to another at the GB interface. When grains 
deform, GBS is an inevitable process as a result of the 
relative movements and rearrangements of the grains 1 . 
GBS is a strongly inhomogeneous mode of deformation 
localized at a very narrow interface layer thus creating 
very strong shear forces. This, together with the fact that 
the shear strength of a general GB between grains of high 
missorientation angle (high-angle GB) is much lower 
than the shear strength of the perfect crystal 2- ’ leads to 
GBS being an efficient deformation mode. During GBS 
the load transfer between the sliding surfaces is 
significantly reduced and the load is redirected to other 
places more resistant to sliding. This sliding creates 
redistribution of the load and appearance of stress 
localization in the microstructure which, in the absence 
of an efficient accommodation mechanism, can lead to 
void formation and microcracking. 

In this paper, a multiscale modelling strategy is 
used to study the effect of GBS for stress localization in 
a polycrystalline microstructure. Although they reveal 
the system behavior at atomic level resolution, MD 
calculations of large systems of grains are 
computationally prohibitive to perform, and a more 
efficient analysis technique is sought. Finite element 
method (FEM) analysis is an obvious choice, however, 
the FEM models cannot a-priori simulate the 
deformation mechanisms found within the systems of 
grains. Thus, the FEM models must be tuned in order 
to reproduce the stress localization observed in the 
atomistic simulations. For this purpose, a 
polycrystalline model of bimodal grain-size distribution 
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was constructed convenient to study by both, MD and 
FEM simulations. The MD simulations revealed the 
behavior of this model with atomic-scale details. 
Additional MD simulations on a bicrystal model were 
performed to extract the elastic and yield properties of 
the GBs presented in the model. These GB properties 
were then used to tune the FEM model to reproduce 
closely the behavior of the test model to match with the 
MD simulation results. The FEM model can then be 
used as part of a multiscale modelling strategy to 
extrapolate the MD-derived information to larger scales 
to study the failure properties of polycrystalline 
materials. 

The paper is constructed as follows. First, details of 
the configuration model used in both, MD and FE 
simulations are discussed. The specific features of the 
MD model are presented next. Further, the yield strength 
of the GB is defined by using MD simulations on a 
bicrystal. This yield strength is used as an input 
parameter to the FEM model. Next, the FEM version of 
the model is discussed and the results of the FEM 
simulation are compared with those from the MD 
simulation. Finally a large-scale FEM simulation of a 
typical microstructure of 100 grains is analyzed to reveal 
the stress distribution due to grain-boundary sliding in a 
more general polycrystalline system. 

The Configuration Model 

To investigate the effect of GBS on the stress 
distribution in a polycrystal of inhomogeneous grain size 
a configuration model, presented in Figures 1(a) and (b) 
was used. This configuration was used in both, MD and 
FEM simulations with the additional goal to tune the 
parameters of FEM to reproduce as close as possible the 
MD simulation. The model is a periodically repeating 
unit of a four-grain octagon-square configuration which, 
when replicated in x- and y-directions represents a 
polycrystal of a bimodal grain-size distribution. The 
advantages of this configuration are the following: (i) the 
model allows exploration of the effect of large variations 
in the ratio, a:b, between the grain sizes (Figure 1(b)), 
which for this study was set to 1:4; (ii) the model is 
both simple enough to allow for MD simulation, and 
informative enough to be used in the FEM simulation to 
compare with the MD results; (iii) the configuration has 
a four-fold symmetry (against 90° rotation) which further 
simplifies the analysis. In both, MD and FEM 
simulations, this configuration is deformed by applying 
stress or strain, respectively, along the x- direct ion. The 
deformation conditions are chosen such that the grains 
deform elastically while the GBs deform plastically by 
GBS. While in the MD model, an interatomic potential 
representative for aluminum (Al) is used, to isolate the 
effect of crystal anisotropy on the stress distribution an 
isotropic media is used for the FE simulation. Though 
this would result in some inevitable discrepancy in the 
results, it would help to identify the amount of stress 
distribution due to GBS alone. Fixed boundary 
conditions, shown in Figure (b), which for an isotropic 


media mimic the periodic boundary conditions in the 
MD simulation, shown in Figure (a), were applied in the 
FE simulation. 
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Figure 1: (a) A snapshot of the thermally equilibrated 
polycrystalline microstructure designed to study the 
effect of GBS on the stress distribution in a polycrystal 
of bimodal grain-size by MD simulation. The GBs are 
shown as darker (blue) lines of atoms which were 
identified to be in disordered (non-crystalline) local 
environment. The dashed lines mark the actual system- 
box dimensions while, for better visualization the 
microstructure is slightly extended using periodic 
boundary conditions, (b) Repeating unit used in FE 
simulations to model grain structure shown in (a). As 
there is no characteristic length-scale in the FE model, 
the dimensions in (b) are relative, characterized by a:b = 
1:4 consistent with (a). 


The Molecular-Dynamics Model 

To investigate the effect of GBS on the stress 
distribution in a polycrystal of inhomogeneous grain- 
size, atomistic MD simulations were used where the 
only predefined input quantity is the interactive potential 
between individual atoms which, following Newtonian 
dynamics, defines their motion. Embedded-atoms- 
method (EAM) many body potential 4 fitted to reproduce 
closely the elastic and thermodynamic properties of a 
perfect A1 crystal was applied. A textured or columnar 
microstructure model with periodic boundary conditions 
in all three dimensions to mimic bulk conditions 5 was 
ideal for this study. While providing a fully three- 
dimensional (3d) treatment of the underlying physics, 
this model makes possible to simulate relatively large 
grains, because only a few lattice planes need to be 
considered in the periodically repeated texture direction. 
The texture z-axis is along [110] crystallographic 
orientation which enables a realistic dislocation-slip 
dynamics along six available slip systems in two {111} 
slip planes in each grain 6 . These slip systems are 
adequate to accommodate any two-dimensional 
deformation in the x-y plane of the simulated structure. 
The grain-boundaries, formed between the grains in this 
model are [110] tilt GBs which have been studied 
extensively in the past 15 years both, experimentally 7 " 10 
and by simulations 9 " 13 . In spite of some drawbacks 
inherent in this model when compared to a fully 3d 
structure, as discussed in reference 14, the model was 
proven remarkably successful in predicting the highly 
unexpected deformation twinning in nanocrystalline 
A1 “ , and the transition from partial to perfect slip in 
nanocrystalline fee metals 617 . Both of these phenomena 
have since gained a solid experimental support " . 

Specifically for this study, a four-grain octagon- 
square configuration was used as a test model. A 
snapshot of the microstructure after thermal equilibration 
at 100 K is given in Figure 1(a). The configuration 
represents a periodically repeating unit of large octagonal 
grains of size, b = 80 nm, encapsulating small square 
grains of size, a = 20 nm. The thickness of the system 
in the texture direction was set equal to 5 (1,1,0) atomic 
planes (1.43 nm), as in the previous columnar 
models 5 ’ 6 ' 1516 . Within these dimensions, the system 
contains 1.7 million atoms. By rotating one of the two 
octagonal grains at 90° relative to the other (the 
crystallographic orientations are given in Figure 1) the 
four 45° inclined diagonal GBs become 90° 299 
symmetric-tilt GBs (STGB) for which the atomic 
structure in A1 is known from the literature 10 . The two 
types of 299 STGBs with (5,5,7) and (7,7,10) GB 
plane, (see Figure 1(a)), are structurally very similar and 
practically undistinguishable from the point of view of 
their mechanical properties, which preserves the four-fold 
symmetry (against 90° rotation) of the structure. In 
addition, the two small square grains are rotated at 45° 
relative to the lar^e grains thus forming high-angle 
asymmetric tilt GBs ’ (Figure 1(a)). 


To initiate GB sliding the microstructure was loaded 
with uniform tensile stress of 1 GPa along the x-axis in 
Figure 1(a). The stress was applied using Parrinello- 
Rahman 20 constant-stress technique combined with Nose- 
Hoover thermostat 21 for constant-temperature simulation. 
This stress is well below the threshold stress of about 2 
GPa 6 ' 14 needed to start nucleation of dislocations from 
the GBs, which are the only possible dislocation sources 
in grains of a nanometer size 14 ' 15 . Thus, the dislocation 
activity, as a possible accommodation mechanism was 
readily suppressed. By running the simulations at very 
low temperature of 100 K (a value still high enough to 
avoid quantum-mechanical effects existing in a real 
structure that cannot be captured by the classical MD 
technique) GB diffusion, as another possible 
accommodation mechanism, was also eliminated. Also, 
low temperature prevents grain-growth, a process which 
otherwise will be very strong in such a system of a 
bimodal grain-size distribution. Restricting the 
dislocation activity and GB diffusion helps to map the 
MD model better to our FEM simulations which, at the 
present stage, do not include these mechanisms. 

During loading the system reached an equilibrium 
strain of 1.23%. The strain was elastic within the grains 
and plastic at the GBs thus initiating GB sliding. At 
this stage, to present the stress distribution in the 
system, two-dimensional stress-maps were created. The 
two-dimensional resolution of the maps is an area of 6x6 
lattice parameters (2.43 x 2.43 x 1.43 nm 3 volume 
including 432 atoms) over which the local stress has 
been averaged. The stress maps were stored every 1 ps 
over a period of 20 ps after reaching equilibrium (approx. 
40,000 MD steps) and then averaged in time to smooth 
out the fluctuations always present in a system of such 
small size. These averaging procedures in space and time 
give a good estimate of the local stress calculated by the 
virial-expansion technique in Parrinello-Rahman stress- 
calculation" 0 . All the stress-maps for the three stress 
components in the x-y plane of deformation (i.e., normal 
to the texture z-axis), o xx , o yy and a xy , show distinctive 
stress distribution and stress localization which were 
then compared with an FEM simulation of the same 
configuration, as will be discussed below. 

Shear strength of a grain-boundary 

Knowledge of the shear strength of the GBs 
presented in the system is crucial for understanding the 
stress distribution due to GBS. This shear strength is 
also required as an input for the FEM model discussed 
below. For this purpose, a separate MD simulation of a 
bicrystal cubic system, presented in Figure 2, was used. 
The same interatomic potential 4 as in the octagon-square 
model was applied. 

The two crystals were crystallographically oriented 
in such a way as to form a (7,7,10) 299 STGB (Figure 
2). This type of GB, together with the very similar 
(5,5,7) 299 STGB, were expected to take most of the 
GBS in the octagon-square model (Figure 1) due to the 



maximum resolved shear stress on the planes inclined at 
45° to the tensile direction. In the simulation, the system 
box with periodic boundary conditions was allowed to 
shear in all three directions in addition to expansion or 
contraction making it possible to estimate the shear 
strain directly. The system was loaded with a symmetric 
shear stress ( Oxy=Oyx) to prevent torque forces and to 
initiate GBS in the direction consistent with the sliding 
direction between the octagons in the octagon-square 
model. In addition, a control simulation was performed 
on a perfect crystal of the same size and orientation 
([5,5,7]/[7,7,l 0]/[ 1,1,0]) to obtain the shear strain of a 
perfect crystal under the same simulation conditions. The 
perfect-crystal strain was then subtracted from the 
bicrystal strain and the resulting strain, due to the GB 
alone, is presented in Figure 3 as a function of load for 
three different temperatures. 

Following Schiotz et al.” the yield stress was 
defined as the stress where the strain starts to depart from 
linearity, a definition convenient for use in MD 
simulations. For the temperature of T = 100 K, used in 
the simulation of the octagon-square model, the obtained 
yield stress is a y = 0.2 GPa. As expected, the yield 
stress decreases with increasing the temperature of the 
system becoming 0.14 GPa at 200 K and 0.12 GPa at 
300 K. For comparison, the theoretical shear strength of 
a perfect A1 crystal was estimated at around 2.8 GPa’ 3 . 
Thus, the GBs are shown to be more than an order of 
magnitude weaker to shear than the perfect crystal grains. 
Applying a tensile load of 1 GPa along the x-axis in the 
octagon-square model (Figure 1(a)) results in a resolved 
shear stress of 0.5 GPa along the diagonal 299 GB 
planes. This stress is more than twice as large as the GB 
yield stress, as defined here (see Figure 3), but still well 
below the shear strength of the perfect crystal. This 
ensures a regime of deformation where the GBs would 
experience plastic sliding, while the crystal grains would 
deform elastically. 

FEM simulation 

For the finite element simulations, a code “FRANC” 
developed at Cornell University 24 was used and 
specifically designed to study crack initiation and 
propagation in metallic polycrystals with explicit 
representations of grains and grain boundaries. Each of 
the grains in the model was given isotropic, elastic 
material properties. The properties used were Young’s 
modulus, E = 74.8 GPa, and Poisson’s ratio, v = 0.346. 
These values were calculated averages’ 3 from the 
anisotropic elastic constants of the interatomic potential 4 
used in the MD approach. To accurately account for the 
entire structure, a repeating unit, similar in geometry to 
the one used in the MD simulation was chosen (Figure 
1(b)). As in the MD model, the ratio of the lengths of 
the grain boundaries is a:b=YA. Periodic boundary 
conditions are necessary to accurately compare results to 
those obtained from MD simulations. The boundary and 
loading conditions shown in Figure 1(b) are consistent 
with periodic boundary conditions for an isotropic 
material 


In this work, cohesive zone models (CZM) were 
used to characterize grain boundary behavior’ . Cohesive 
zone models assume cohesive interactions of the material 
around a grain boundary and permit the appearance of 
fracture surfaces in a continuum 27 . To simulate the 
sliding behavior along grain boundaries, two 
independent cohesive zone models are chosen for the 
normal and shear components of the traction and 
displacement. To permit sliding, a perfectly plastic 
relationship is chosen for the shear model with yield 
stress of 0.2 GPa; this stress is consistent with the 
numbers obtained from the MD criterion discussed 
previously. To restrict opening, the normal CZM is 
specified as having linear elastic behavior with high 
stiffness. The analysis is performed for 1.23% applied 
strain. This applied strain corresponds to the 1 GPa 
applied stress in the MD simulations. 



Figure 2: Bicrystal cubic microstructure used to get the 
shear strength of 299 STGB formed in the octagon- 
square system shown in Figure 1. 



Figure 3: GB shear strain vs stress extracted from the 
bicrystal model shown in Figure 2. The straight lines, 
tangential to the curves at their origin, help to mark the 
linear elastic regime. The arrows show the yield stress 
for each of the three simulated temperatures defined 2 ’ as 
the stress where the strain starts to depart from linearity. 
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Figure 4: (a-c) Stress contours from the FEM model and the corresponding stress maps (d-e) from the MD model for o xx , 
o yy , and o xy stress components. Positive and negative stress are indicated in red and blue, respectively, relative to the 
average background stress defined as white. In (a), (b), (d) and (e) positive stress is defined as the stress in compression and 
negative as the stress in tension compared to the background stress. In (c) and (f) positive and negative shear corresponds to 
shear directions at the diagonal GBs defined by Figure 5. The stress localization of o yy at the GBs of the small grain, 
marked with (1) in (b) and (e) are significant. Marked with (2) in (c) and (f) are also the shear stress concentrations created 
at the triple junctions. A few GB dislocations that appeared in the MD simulation and which create disturbances in the GB 
stress field are indicated with (3) in (f). 


Results 

Applying a uniform tensile stress of 1 GPa creates 
resolved shear stresses larger than the GB yield stress, 
but lower than the yield stress of the grains, producing 
plastic sliding at the GBs and elastic deformation in the 
grain interiors. GBS reduces the load transfer from one 
grain to another and redistributes the stress. The a xx , o yy , 
and o X y stress distributions for the octagon-square model 
are shown in Figures 4(a-f). The results from the FE 
simulations Figures 4(a-c) are compared with the stress 
maps obtained from the MD simulation Figures 4(d-f). 
Both models show that the small square grains take 
most of the load produced from the tensile stress (see 
Figures 4 (a) and (d)) and from the Poison contraction 
perpendicular to the tensile x-direction (see Figures 4 (b) 
and (e)). This contraction creates a strong tensile stress at 


the sides of the square grains normal to the tensile 
direction (regions marked (1) in Figures 4 (b) and (e)). 
The close similarity between MD and FEM stress 
distribution and the fact that the FEM model was 
isotropic indicates that this stress distribution is solely 
due to GBS. Through GBS the grains become decoupled 
in their elastic deformation. At the GBs, where only two 
grains meet, this decoupling creates sliding, but at the 
triple junctions where three grains meet strong 
incompatibility stresses are produced. This is seen in 
both FEM and MD results (Figures 4 (c) and (f». 

The most direct evidence for the presence of GBS is 
seen on the MD o xy -stress map in Figure 4 (f). The four 
diagonal GBs showed strong shear stress in the two 
opposite directions (with stress values presented in red 
and blue). Figure 5 gives a schematic representation of 
the grain sliding responsible for the observed shear stress 





distribution. Grains (1) and (3) are displacing towards 
each other compressing the small grain (5) in the y- 
direction, while grains (2) and (4) are displacing apart 
and stretching grain (5) in the x-direction. The 
displacement of the large grains results in sliding at the 
GBs, schematically shown in Figure 5 as couples of 
oppositely directed arrows on both sides of the diagonal 
GBs. The direction of each arrow corresponds to the 
local displacement of each grain at the GB. This 
displacement produces strong local shear stress at the 
diagonal GBs as seen in Figure 4(c,f). 

In spite of the overall good match between the stress 
distribution from MD and FEM simulations, there are 
also significant mismatches that can be discovered after a 
careful examination of Figures 4(a-f). For example, in 
the MD simulation there is no full symmetry of the 
stresses produced on the two small grains, and there are 
some deviations in the stresses on the large grains too 
(most pronounced in Figure 4(d)). Possibly, the main 
reason for these deviations is the small, but noticeable 
anisotropy in the MD model which results from the 
interatomic potential. This potential was fitted 4 to 
reproduce closely the elastic constants, including the 
anisotropy, of monocrystalline Al. An anisotropic 
FEM simulation with the proper periodic boundary 
conditions is needed to get better comparison. 

Another source of mismatch between the two 
models is the inevitable existence of GB dislocations in 
the MD model. In a polycrystalline MD model, 
absolutely perfect GBs are impossible to sustain during 
loading, even if they were produced initially. The same 
is true for a real material, but at this stage GB 
dislocations are not considered in the FEM model. The 
GB dislocations have long-range elastic fields that can 
significantly modify the stress distribution in close 
proximity to the GBs thus altering it from the FEM 
result. Such an effect produced by GB dislocations is 
marked as (3) in Figure 4(f). Also, the phonon waves, 
always existing in finite temperature MD simulation, 
create periodic modulations in the MD stress maps 
Figures 4 (d-f) which are non-existent in the FEM stress 
distributions (Figures 4 (a-c)). Even averaging of the 
stress over 20 ps time interval was not enough to 
smooth these fluctuations enough. One reason may be 
that the periodic boundary conditions imposed on the 
system, together with its high configuration symmetry, 
create stationary waves which cannot be smoothed out 
over such a short period of time. All these factors make 
the MD maps shown in Figure 4 much more rich and 
complicated compared to the substantially more idealized 
FEM results. 

The close quantitative match between the stress 
distribution obtained from MD and FEM simulations of 
the octagon-square configuration (Figure 4) gives 
confidence that implementing the properly parameterized 
CZM elements in the FEM model correctly reproduces 
GBS for a polycrystalline specimen. With this 
confidence, FEM can now be used to simulate a typical 
polycrystalline structure. For this purpose, a sample with 


100 randomly shaped and sized grains, represented by 
Voronoi polygons' 8 , was constructed. The grain 
configuration, along with the boundary and loading 
conditions, are shown in Figure 6. 



Figure 5. A schematic representation of the sliding 
directions of the four octagonal grains (bold arrows) and 
the inserted stresses on the small square grain in the 
middle (dashed arrows). The local displacements of the 
material at the GBs, producing GBS are presented by 
couples of opposite pointing arrows on both sides of the 
boundaries. 
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Figure 6: Material model with 100 randomly shaped and 
sized grains represented by Voronoi polygons 28 . 



GPa 


Figure 7: Stress contour o xx for 100-grain material subjected to 1.23% strain, as shown in Figure 6, with some similarities 
to the idealized model highlighted. 


For 1.23% applied strain, the o xx stress distribution 
is shown in Figure 7. The similarities in grain boundary 
behavior between the randomly configured model and the 
idealized octagon-square configuration are highlighted. 
This larger-scale simulation shows that stress 
localization due to GBS are common within the 
polycrystals and local stress can exceed several times the 
background stress, thus marking spots for a possible 
crack nucleation. It is not surprising -0 that these stress 
concentrations are mostly localized at the triple junctions 
where three GBs meet. 

Summarizing the results from this section, it was 
shown that simulating the behavior of geometrically 
identical polycrystalline structures by MD and FEM 
techniques can serve two purposes. First, it helps to tune 
the parameters of the larger-scale FEM simulation to the 
results from the lower-scale MD simulation, thus 
incorporating the effects of atomistic mechanisms, 
presented in MD, into a continuum model. Second, such 
a parallel use of two very different techniques can help to 
identify more unambiguously the underlying mechanism 
behind a certain phenomenon. In the presented study, 


reproducing a stress distribution by isotropic FE 
simulation similar to the one obtained from MD 
simulation where anisotropy is inherent from the 
interatomic potential and cannot be turned off, proves 
that GBS played the major role in this case. As a result, 
it was shown that even in a material of weak anisotropy, 
such as Al, substantial stress localization can be 
produced due to GBS. 

Concluding Remarks 

In the presented study a combined implementation 
of two modelling techniques, MD and FEM, has been 
demonstrated working at different length-scale 
resolutions on two geometrically equivalent systems. 
The results, after comparing these two simulations, are 
very instructive. First, a method was successfully 
applied to match the parameters in the constitutive 
equations (in this case the CZM) used in the FEM 
model with quantities extracted from the MD simulation 
(like GB yield stress). By this method the FEM was 


tuned to the MD simulation. Second, in an MD 
simulation where the resolution is at the level of 
individual atoms, it is usually very difficult to identify 
unambiguously the underlying mechanism responsible 
for a given process. By reproducing the MD simulation 
experiment by FEM model, where the mechanisms are 
implemented through constitutive equations and can be 
easily controlled, the underlying mechanism can be 
revealed with higher certainty. For example, in the 
octagon-square model, presented here, the stress 
distribution obtained by MD simulation cannot be 
attributed readily to GBS, as there is also anisotropy 
which can produce analogous effects. After reproducing 
the similar simulation experiment by FEM using 
isotropic material parameters, but allowing for GBS with 
the equivalent sliding resistance as in the MD model and 
comparing the results it becomes possible to identify 
what part of the stress distribution is due to GBS alone, 
and what part is due to other sources like anisotropy and 
GB dislocations. 

In this study it has been also demonstrated that 
now, because of the advanced computational technology, 
high-end MD simulations can overlap well with FEM 
simulations on yet small, but still informative enough 
systems to allow direct comparison of the results. This 
possibility suggests a new way of bridging length and 
time scales, which is an alternative to the direct multi- 
scale models, where the simulated system is divided into 
regions of different length- and time-scales which are 
treated with different simulation techniques 
simultaneously 3 '’. Instead of trying to combine the two, 
inherently very different simulation techniques, into one 
multi-scale modelling program, it is now possible to use 
the overlapping region to tune the larger scale model to 
the lower-scale one by direct comparison. This approach 
has the advantage that it preserves the uniqueness and the 
integrity of each of the models and avoids creating 
overcomplicated software codes with very restricted 
applicability. This approach also avoids the artifact of 
the interface between regions of different length- and 
time-scales existing in the direct multi-scale models 30 . 

Finally, the large-scale implementation of the MD- 
tuned FEM simulation on a typical 100 grain 
microstructure showed that GBS, by redistributing the 
load in the system, can produce stress localization 
several times higher than the background stress under 
external uniform load. These stress localization, under 
suitable conditions, may become nuclei for 
microcracking and void formation. 
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